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Abstract 

The  present  paper  describes  numerical  modelling  of  the  radiative  heat  transfer  process  in  the  module  chamber  of  an  internal  indirect  reforming- 
type  SOFC.  The  ability  to  do  internal  reforming  is  one  of  the  characteristics  of  high-temperature  fuel  cells,  SOFC.  As  in  any  high-temperature 
system,  radiative  heat  transfer  is  important.  In  this  article,  heat  transfer  between  the  fuel  reformer  surface  and  all  other  surfaces  facing  the  reformer 
surfaces  is  modelled.  Governing  equations  for  radiative  heat  transfer  are  described  using  Hotter  s  zone  method.  The  resulting  radiation-conduction 
conjugate  heat  transfer  problems  are  numerically  solved  with  a  combination  of  Gauss-Seidel  and  Newton-Raphson  methods.  The  steam  reforming 
reaction  occurring  inside  the  fuel  reformer  is  described  using  Achenbach  model.  The  obtained  results  indicate  that,  for  the  development  of  effective 
indirect  internal  reforming,  the  position  of  the  reformer  in  the  module  chamber  and  emissivity  of  the  surfaces  of  the  reformer,  cell  and  other 
elements  in  the  SOFC  module  all  play  a  key  role. 

©  2008  Elsevier  B.Y.  All  rights  reserved. 
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1.  Introduction 

Due  to  their  high  efficiency  and  the  low  level  of  pollution 
they  emit  into  the  environment,  solid  oxide  fuel  cells  (SOFCs) 
have  the  potential  to  become  one  of  the  most  important  energy- 
conversion  devices.  A  solid  oxide  fuel  cell  consists  of  two 
porous  ceramic  electrodes  (the  cathode  and  anode)  separated 
by  a  solid  ceramic  electrolyte.  A  typical  SOFC  uses  solid  oxides 
such  as  yttria-stabilized  zirconia  (YSZ)  for  the  electrolyte  [1]. 
It  can  use  hydrogen,  carbon  monoxide,  or  hydrocarbons  as  fuel 
and  air  or  oxygen  as  oxidants.  The  operation  temperature  of 
an  SOFC  is  700-1 100  °C  and  operates  at  atmospheric  or  ele¬ 
vated  pressures  [2].  Hydrocarbons  are  the  most  suitable  fuels 
for  SOFCs  to  achieve  high  efficiency  from  natural  resources 
[3].  For  hydrocarbon-based  fuel,  three  types  of  fuel  conversion 
can  be  considered  in  reforming  reactions:  an  external  reform¬ 
ing  system,  an  indirect  internal  reforming  system  and  a  direct 
internal  reforming  system.  High-temperature  SOFC  eliminates 
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the  need  for  an  expensive  external  reforming  system.  The  pos¬ 
sibility  of  using  internal  reforming  is  one  of  the  characteristics 
of  high-temperature  fuel  cells,  SOFC.  Strong  endothermic  fuel 
reforming  reactions  can  be  thermally  supported  by  the  transfer 
of  the  heat  generated  due  to  the  sluggishness  of  electrochemical 
reactions,  diffusion  of  participating  chemical  species  and  ionic 
and  electric  resistance.  However,  for  high-temperature  opera¬ 
tion,  thermal  management  of  the  SOFC  system  becomes  an 
important  issue.  To  carry  out  thermal  management  properly, 
detailed  modelling  and  detailing  of  numerical  analyses  of  the 
phenomena  occurring  inside  the  SOFC  system  are  required. 

Thermal  models  of  fuel  cells  that  use  heat  transfer  correla¬ 
tions  for  assumed  fluid  flow  configurations  have  been  designed 
in  the  past  [4].  The  Suzuki  group  [5,6]  proposed  a  quasi-two- 
dimensional  model  for  numerical  studies  on  the  performance  of 
a  single  tubular  SOFC  cell  under  practical  operating  conditions. 
Their  model  takes  into  account  the  air  and  fuel  flow  velocity 
fields,  ohmic  and  thermodynamic  heat  generation,  convective 
heat  transfer,  mass  transfer  of  participating  chemical  species, 
including  the  electrochemical  processes  and  the  electric  poten¬ 
tial  and  electric  current  in  the  electrodes  and  electrolyte  [5].  The 
results  from  this  model  were  applied  to  the  performance  anal¬ 
ysis  of  a  tubular  SOFC/micro-gas  turbine  hybrid  system  [6]. 
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Nomenclature 

A  area  (m2) 

Cp  specific  heat  ( J  kg- 1  K- 1 ) 

d  thickness  (m) 

E  black-body  flux  density  (W  m-2) 

Hj  incident  radiant  flux  density  (W  m-2) 

lij  distance  between  ith  and  jth  elements  (m) 

m  molar  fraction 

Mj  molar  mass  of  j  component  (kg  mol-1) 

h  molar  flow  rate  (mol  s - 1 ) 

P  pressure  (Pa) 

q  heat  flux  (W) 

R  reaction  rate  (mol  m-  3  s  - 1 ) 

Ru  universal  gas  constants,  8.314472  (J mol-1  K-1) 
SiSj  direct  exchange  area  (m2) 

SiSj  total  interchange  area  (m2) 

Sr  heat  generation  coming  from  reforming  reaction 

(Wm-3) 

T  temperature  (K) 

U  velocity  (ms-1) 

W  radiosity  (W  m-2) 

x,  y,  z  axial  coordinates 

Greek  symbols 

8  Kronecker  delta 

s  emissivity 

£o  porosity 

0  angle 

k  wavelength  (m) 

kj  thermal  conductivity  (W  m-1  K-1) 

/ 1  dynamic  viscosity  (Pa  s) 

p  reflectance 

po  density  (kg  m- 3 ) 

o  Stefan-Boltzmann  constant,  5.6667  x  10-8 

(Wm-2  K-4) 

(p  view  factor 

Subscripts 

B  bottom  of  geometry  system 

eff  effective 

f  fluid 

i,  j  ordinals  numbers 

max  maximum  value 

o  ambient 

R  reformer 

sd  solid 

sh  shift  reaction 

st  steam  reforming  reaction 

S  side  walls 

T  top  of  geometry  system 


A  steady- state  model  of  an  indirect  internal  reforming  solid 
oxide  fuel  cell  was  developed  by  Aguiar  et  al.  [2]  to  anal¬ 
yse  the  thermal  coupling  of  the  reforming  and  electrochemical 
reactions.  In  their  analysis,  the  configuration  was  an  annular 


design  where  the  fuel  call  was  constructed  around  a  tubular 
reformer.  Their  model  was  based  on  a  heterogeneous,  two- 
dimensional  steam  reforming  packed-bed  reactor  model  coupled 
with  a  one-dimensional  SOFC  one.  Temperature  distributions  in 
three-dimensional  models  of  planar  SOFC  have  been  achieved 
by  coupling  the  electrochemical  effects  with  the  solution  of 
the  Navier- Stokes  equations  of  motion  and  energy/mass  con¬ 
servation  [7] .  As  in  any  high-temperature  system,  radiative  heat 
transfer  is  important.  However,  most  existing  models  neglect 
its  effects  [8].  Yakabe  et  al.  [9]  found  that  including  radiative 
heat  exchange  within  the  flow  channels  resulted  in  a  flatter  dis¬ 
tribution  of  the  temperature  profile  along  the  fuel  cell.  Murthy 
and  Fedorov  [8]  demonstrated  radiation  heat  transfer  effects  for 
accurate  prediction  of  the  temperature  field  and  fuel  cell  out¬ 
put  voltage.  Damm  and  Fedorov  [10]  have  recently  reviewed 
radiative  heat  transfer  in  SOFC  components  and  materials. 

Experiments  have  been  carried  out  by  Achenbach  and  Rien- 
sche  [11]  to  determine  the  kinetics  of  the  methane/steam 
reforming  process  in/on  anode  materials  of  a  solid  oxide  fuel 
cell.  Their  tests  [11]  showed  that  H2O  partial  pressure  has  no 
catalytic  effect  on  the  reaction.  They  have  also  proposed  an  equa¬ 
tion  for  the  reforming  kinetics  in  the  form  of  the  Arrhenius-type 
independence  of  the  H2O  partial  pressure  and  proportional  to 
the  methane  (CH4)  partial  pressure  [11].  Their  methane/steam 
reforming  kinetics  model  [11]  has  been  adopted  to  numerical 
analysis  of  a  planar  solid  fuel  cell  stack  [12]  as  well  as  to  a 
tubular  solid  oxide  fuel  cell  with  accompanying  indirect  inter¬ 
nal  fuel  reforming  [6].  An  alternative  model  of  the  methane 
reforming  reaction  rate  was  proposed  by  Odegard  et  al.  [13]. 
This  model  [13]  was  used  in  numerical  analyses  of  the  char¬ 
acteristics  of  tubular  SOFC  with  internal  reformer  systems  by 
different  research  groups  [14,15]. 

The  present  paper  describes  a  three-dimensional  numerical 
model  of  heat  transfer  process  in  the  module  chamber  of  an 
internal  indirect  reforming-type  SOFC.  Heat  transfer  between 
the  fuel  reformer  surface  and  all  other  surfaces  facing  the 
reformer  surface  is  modelled.  The  Achenbach  model  [12]  is  use 
to  describe  the  steam  reforming  process  occurring  inside  the  fuel 
reformer  and  resulting  radiation-conduction-convection  con¬ 
jugate  heat  transfer  problems  are  numerically  solved  by  the 
Gauss-Seidel  and  Newton-Raphson  methods.  This  model  can 
be  applied  for  different  configurations  and/or  characteristics  of 
the  system  components. 

2.  The  formulation  and  numerical  scheme 

The  geometry  of  the  problem  to  be  considered  is  shown 
schematically  in  Fig.  la  (a  detailed  description  on  the  system 
discussed  in  this  paper  may  be  found  in  Ref.  [27]).  The  surfaces 
in  the  system  are  assumed  to  be  grey  and  the  gases  in  the  SOFC 
chamber  are  assumed  to  be  transparent.  The  supplied  fuel  to  the 
system  is  the  methane  (CH4).  The  gaseous  fluids  are  assumed 
to  be  Newtonian  and  the  flow  of  fuel  inside  the  fuel  reformer  is 
assumed  to  be  laminar,  steady  and  occurring  in  one  direction. 
The  effect  of  the  buoyant  convection  is  not  considered.  All  chem¬ 
ical  components  of  working  fluids  in  the  system  are  considered 
to  be  the  ideal  gases  and  the  temperature  of  in-flow  fuel  is  taken 
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Fig.  1.  Configuration  of  the  SOFC  system  (A)  and  the  grid  used  for  numerical  computation  (B). 


to  be  instantly  approaching  the  temperature  of  the  reformer.  The 
analysis  of  heat  transfer  between  the  reformer  and  SOFC  cham¬ 
ber  illustrates  that  more  than  98%  of  energy  is  transported  by 
radiation,  less  than  1 .9%  of  energy  is  transported  by  convection 
and  less  than  0.02%  of  energy  is  transported  by  conduction. 


2.7.  Radiative  heat  transfer  model 


The  geometry  of  the  radiative  heat  transfer  model  to  be 
considered  is  shown  schematically  in  Fig.  lb.  The  govern¬ 
ing  equations  to  be  solved  are  energy  balances.  Assuming  the 
absence  of  free  and  forced  convection  inside  the  closed  surface 
system,  the  rigorous  formulation  of  energy  balances  on  various 
surfaces  yields  an  integral  equation  which  can  be  written  as 


Wit/cos  Oij  cos  OjiSxj 


nrf 


dAjdX 


SkiEudX  + 


(1) 


where  d A/  is  the  area  of  the  free  surface  element,  Oq  is  the  angle 
between  rq  and  the  normal  to  d A/,  and  q  are  the  leaving 

flux  density  (radiosity),  wavelength,  black-body  flux  density, 
emissivity  and  heat  flux,  respectively. 

The  energy  balance  on  the  surface  A/  (Eq.  (1))  contains  var¬ 
ious  terms:  absorption  of  incidental  radiation  from  all  other 
surface  elements,  emission  by  surface  A/  and  the  heat  flux 
extracted  through  surface  A/.  Considering  an  enclosure  of  any 
shape  with  walls  varying  in  temperature  in  any  manner  and 
assuming  that  the  walls  are  grey,  Eq.  (1),  based  on  the  Hottel 
zone  method  [16],  is  written  as 


SjSiEj  —  AiSiEi  +  Atqt  =0,  i  =  1,  . . . ,  m  (2) 

7=1 

where  Ei-ol 4,  Ei  is  the  black-body  flux  density,  T 
is  temperature,  a  is  the  Stefan-Boltzmann  constant 
(cr  =  5.6667  x  10-8Wm-2  K-4),  Si  is  the  emissivity  of 


the  ith  surface,  and  SjSi  is  the  total  interchange  area  between 
surface  zones  st  and  sj. 

The  radiosity  Wj  can  be  expressed  as 

wj  =  sjej + PjHj 

where  Hj  and  pj  are  incident  radiant  flux  density  and  reflectance, 
respectively. 

The  value  of  the  total  interchange  area  SjSi  between  surface 
zone  Si  and  surface  zone  sj  is  determined  by 


SjSi  =  (  )  G iWj  -  SijSi) 


(3) 


where 


SiWj  = 


AjSj 

Pi 


D 


D  represents  the  determinant  of  the  matrix  [sis]  — 
8ij(Ai/pi )],  8q  is  the  Kronecker  delta  and  D q  is  the  cofactor 
of  matrix  [sis]  —  8q{Ai/  pi)]. 

The  direct  interchange  area  s]s]  or  the  view  factor  (pq,  between 
the  surface  zones  Si  and  sj  may  be  evaluated  by  integration  of 
the  equation 


[(cos  Oij  cos  0ji)/(7Trfj)]dAidAj 


(4) 


where  sis]  =  Ai(pq,  cpq  is  the  view  factor,  and  A/  the  area  of  the 
ith  surface. 

For  a  closed  system  of  surfaces  containing  transparent  gas, 
the  following  relation  should  be  satisfied: 


sts  j  =  s  jSi  or  Ajcpij  =  A  j<p  /7 
and 

YA°j  =  Ai 

j 


ESiSj  —  £‘A‘ 

j 
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2.2.  Mathematical  model  of  internal  indirect  reforming 
SOFC  system 


2.2.1.  Heat  transfer  model  inside  the  reformer 

The  fuel  is  continuously  supplied  to  a  reformer  to  maintain 
appropriate  operating  conditions  for  the  fuel  cells  stack.  The 
most  basic  governing  equations  to  be  solved  by  numerical  mod¬ 
elling  are  for  the  velocity,  temperature  and  concentration  fields 
in  the  internal  reformer.  They  are  different  among  two  kinds  of 
areas:  solid  area  and  porous  area. 

For  the  solid  part  (the  wall  of  the  fuel  reformer)  the  gov¬ 
erning  equation  to  be  solved  is  considered  to  be  the  following 
three-dimensional  heat  conduction  equation  without  internal 
heat  generation: 


dT  (  dT\ 
dx  \  W  3 x  ) 


=  0 


(5) 


where  A,w  [W  m-1  K-1]  is  the  conductivity  of  the  solid  wall  of 
the  reformer. 

The  heat  flux  in  Eq.  (2)  can  be  expressed  as 


<li  = 


Tj  -  Tj-i 

Ri 


(6) 


where  Rj  =  dfkw,  Ti  is  the  temperature  of  the  ith  surface  and  d 
is  the  thickness  of  the  reformer  wall. 

The  microstructure  of  porous  media  is  not  to  be  considered 
directly  in  the  present  model  therefore  for  the  porous  area,  the 
governing  equations  derived  by  the  volume- averaging  method 
are  applied.  In  the  adopted  method  physical  values  are  locally 
averaged  for  a  representative  elementary  volume  [17].  The  rep¬ 
resentative  volume  is  sufficiently  larger  than  the  scale  of  the  fine 
structure  of  the  porous  medium  and  is  sufficiently  smaller  than 
the  scale  of  the  porous  body  itself  [15].  The  heat  transfer  pro¬ 
cess  inside  the  reformer  porous  media  is  described  by  continuity, 
momentum,  energy  and  mass  transfer  equations.  Consequently 
the  following  transport  equations  of  the  averaged  physical  values 
for  laminar  flows  apply 

dppUx  dppUy  dppUz  _ 

dx  dy  dz 


Po  (  dUx  dUx 
Ux  dx  +  Uy  dy 


c0 


+  UZ 


dUx 

dz 


dP  1 

dx  £o 


3 

dx 


li- 


3 Ux 
dx 


3 

+  dz 


11- 


ZUx 

dz 


~  j;Ux  -  j=Uxy/Ut  +  Uj  +  U\  (8a) 


9 

+  dy 

pof , 


p- 


dUx 

dy 


^1  v, 


SUy 

dx 


+  Uy 


3  P  1 

“a - ' - 

dy  £0 


+ 


dz 


fl- 


dUy 

dz 


dUy 

dy 


+  uz 


dUy 

dz 


dl^- 


dUy 

dx 


3 

+ 


P 


K 


Uv 


po£ 

Jk 


p- 


dUy 

dy 


u. 


xfuf+uf+vi  <8b) 


£0 


3 Uz 
dx 


+  U, 


3  P  1 

—  — I - 

dz  £0 


;  fy 

3 

dx 


+  UZ 


dUz 

dz 


/X 


dUz 

dx 


3 

+  dz 


ii- 


dUz 

dz 


“  kUz  _  +  uy +  (8c) 


9 

+  dy 

A)/, 


P- 


dUz 

dy 


,  dT  dT  dT 
p0CP  (  ux-  +  Uy-  +  Uz- 


+  Q  (9) 


,  dP0Yj  .  tT  dpoYj 
x  dx  +Uy  dy 


+  UZ 


dpoYj 

dz 


9 

dx 


+  S1 


(10) 


The  physical  values  in  the  above  equations  represent  the  local 
phase  average  of  the  gas  control  volume,  Ux,  Uy  and  Uz  are 
the  gas  phase  average  x-,  y-  and  z-components  of  the  local  gas 
velocity.  T  is  the  local  average  temperature  over  both  the  gas 
and  solid  phases.  The  £o  =  0.9  and  ^=1.0x  10-7m2  are  the 
porosity  and  permeability  of  the  porous  medium,  respectively. 
The  /=  0.0088  is  the  inertia  coefficient  which  depends  on  the 
Reynolds  number  and  the  microstructure  of  the  porous  medium 
[18].  Aeff  =  {fi(Af  +  (l  —  £o)^s}  [Wm-1K_1]  is  the  effective 
thermal  conductivity,  where  kf  represents  fluid  thermal  con¬ 
ductivity  and  ks  =  lO.OWrn-1  K-1  is  the  solid  phase  thermal 
conductivity  [15].  D/W>eff  =  (1  —  Vl  —  £o)Djm  [m2  s_1]  is  the 
effective  mass  diffusivity  of  species  j.  Yj  is  the  mass  fraction  of 
chemical  species  j  and  Djm  is  the  mass  diffusivity  of  species 
j  in  the  multi-component  mixture  of  gases,  respectively.  All 
gas  components  are  treated  as  ideal  gases  in  the  calculation  of 
density.  In  the  calculations  of  the  thermal  field,  temperature  is 
adapted  as  a  variable  to  solve  but  the  thermodynamic  properties 
(specific  heat  at  constant  pressure,  enthalpy,  Gibb’s  free  energy) 
and  transport  properties  (viscosity,  heat  conductivity)  are  treated 
as  functions  of  temperature.  Their  local  values  are  calculated 
at  local  temperature  by  taking  data  from  the  liquids  and  gases 
properties  tables  [19].  Properties  of  mixtures  are  evaluated  with 
mixing  laws  [20]. 


2.2.2.  Mathematical  model  for  fuel  reforming  process 

The  methane/steam  reforming  process  is  widely  known  as 
a  conventional  process  for  producing  hydrogen  [21,22]  and  it 
proceeds  on  catalysts  such  as  nickel-alumina.  In  the  methane 
with  steam  reforming  process,  the  dominant  reactions  are  the 
following  two  [14]: 
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•  fuel  reforming  reaction 
CH4  +  H20  =  3H2  +  CO, 

A  //|300  =  227.5  kJmol- 1  (1 1) 

•  shift  reaction 

CO  +  H20  =  H2  +  co2, 

A^?300  =  — 31.8  kJ  mol-1  (12) 

The  fuel  reforming  reaction  described  by  Eq.  (11)  is  a  slow 
and  highly  endothermic  reaction;  therefore  a  supply  of  thermal 
energy  is  needed  for  the  reforming  reaction.  Also,  how  methane 
is  reformed  inside  the  reformer  is  dependent  on  the  local  con¬ 
ditions,  including  not  only  temperature  but  also  partial  pressure 
of  each  chemical  species  as  well  as  the  density  of  the  cata¬ 
lysts  or  catalysts’  active  area  [15].  The  water-gas  shift  reaction 
described  by  Eq.  (12)  is  a  fast  and  weak  exothermic  reaction  and 
can  be  assumed  to  be  in  equilibrium  at  the  reforming  temperature 
[14]. 

The  reaction  rates  of  the  two  above  reactions  (Eqs.  (11)  and 
(12))  Rst  and  Rs h  are  locally  calculated  as  follows: 

~^i)  (13) 

RT  J 

Rsh  =  K+PCOPR20  -  KshP^lPC02  (14) 


reaction  (Eq.  (12))  is  also  calculated  as  follows: 

5h2  =  3Rst^H2  +  Rs\\Mn2 

(16a) 

Sc  0  =  ^stAfco  RshAfco 

(16b) 

5co2  =  tfshMCo2 

(16c) 

ScH4  =  -^stAfcH4 

(16d) 

5h2o  =  -Rst^u2o  RshMn2o 

(16e) 

The  value  of  mass  production  or  consumption  rate  for  each 
chemical  species  is  introduced  into  the  species  mass  transfer 
(Eq.  (10))  as  a  part  of  its  source  term. 

The  thermodynamic  heat  generation  rate  by  the  reforming 
reactions  (Eqs.  (11)  and  (12))  is  calculated  based  on  the  reaction 
rates,  as  follows: 

2  st  —  A  HS[Rsi 

(17a) 

2sh  =  —  A  HshRsh 

(17b) 

where  A Hst  and  AH sh  are  the  enthalpy  change  accompanied 
with  each  reaction. 

The  heat  conduction  and  convection  equations  were  solved 
numerically  by  finite  volume  method  [24].  The  resulting 
radiation-conduction-convection  conjugate  heat  transfer  prob¬ 
lems  were  subsequently  numerically  solved  by  Gauss-Seidel 
[25]  and  Newton-Raphson  methods  [16]. 


Rst  =  ^cat^CH4PCH4  exp 


Eq.  (13)  is  based  on  the  Achenbach  model  [12]  where 
A cat  is  the  catalysts  active  area  of  steam  reforming  reaction 
in  a  control  volume  of  reformer  [m2m-3],  T  [K]  is  the  fuel 
conversion  temperature,  Ru  =  8.314472  J  mol-1  K-1  is  the  uni¬ 
versal  gas  constant,  Equ4  =  8.2  x  104  J  mol-1  is  the  activation 
energy  of  the  fuel  reforming  reaction  [11],  kcu4  =  4274  x 
10-5  molm-2  Pa-1  s-1  is  the  pre-exponential  factor  [11]  and 
/?ch4  is  the  partial  pressure  of  methane  [Pa].  The  active  area  Acat 
is  an  important  parameter  in  a  sense  that  control  of  its  distribu¬ 
tion  can  be  used  as  a  means  to  change  the  distribution  pattern  of 
the  reformer  temperature. 

and  K~h  denote  the  rate  constants  of  forward  and  back¬ 
ward  water-gas  shift  reactions  and  the  value  of  Rs h  is  determined 
by  a  method  described  by  Lehnert  et  al.  [23].  The  water  shift 
reaction  determined  by  Eq.  (12)  reaches  equilibrium  rapidly; 
therefore  C02,  H2,  CO  and  H20  have  to  satisfy  the  equi¬ 
librium  equation.  Chemical  equilibrium  is  represented  by  the 
equilibrium  constant,  which  is  a  function  of  temperature  and  is 
equal  to  the  ratio  between  the  reactants’  and  products’  partial 
pressures. 


^sh  = 


*1 


K, 


sh 


Pco2PH2  =exp  (_^Pk) 
PcoPu2o  V  RT  J 


(15) 


where  AG°h  is  the  change  of  standard  Gibbs  free  energy  of  shift 
reaction. 

This  equilibrium  constant  given  by  Eq.  (15)  is  introduced  into 
Eq.  (14)  to  calculate  the  rate  of  the  shift  reaction. 

The  mass  production  or  consumption  rate  of  each  chemi¬ 
cal  species  by  the  fuel  reforming  reaction  (Eq.  (11))  and  shift 


3.  Numerical  results 

Numerical  simulation  is  a  potential  tool  for  investigating  the 
flow,  thermal  and  chemical  species  concentration  fields  in  the 
module  chamber  of  an  internal  indirect  reforming-type  SOFC 
process.  In  practice  this  process  involves  radiation  between 
surfaces,  heat  conduction  in  the  solid  area,  forced  convection 
in  the  porous  one  and  the  methane/steam  reforming  process 
proceeds  on  the  catalysts.  However,  applying  numerical  sim¬ 
ulation  enables  the  artificial  exclusion  of  the  forced  convection 
effect  in  the  porous  area  or  the  methane/steam  reforming  pro¬ 
cess  effects  in  order  to  investigate  transport  phenomena  in 
detail.  Therefore,  the  calculations  are  presented  for  two  differ¬ 
ent  problems:  the  radiation-conduction  conjugate  heat  transfer 
problem  without  the  methane/steam  reforming  process  and  the 
radiation-conduction-convection  conjugate  heat  transfer  prob¬ 
lem  with  the  methane/steam  reforming  process.  Numerical 
computations  for  both  problems  were  performed  for  two  differ¬ 
ent  thermal  boundary  tests  (Tests  I  and  II).  Test  I:  the  first  kind  of 
thermal  boundary  was  considered,  i.e.  the  temperature  distribu¬ 
tion  of  the  SOFC  stack  and  the  temperature  of  other  elements  in 
the  SOFC  module  were  assumed  (see  Fig.  1).  Test  II:  the  second 
kind  of  thermal  boundary  was  considered,  i.e.  constant  heat  flux 
from  the  SOFC  stack  and  constant  ambient  temperature  were 
assumed.  Numerical  calculations  were  carried  out  for  different 
positions  of  the  reformer  in  the  module  chamber  and  emissivity 
of  the  surfaces  (i.e.  different  emissivity  of  the  reformer,  cell  and 
other  elements  in  the  SOFC  module).  Table  1  summarized  the 
thermal  boundary  conditions  applied  in  computations. 


Table  1 

Thermal  boundary  condition 


G.  Brus,  J.S.  Szmyd  /  Journal  of  Power  Sources  181  (2008)  8-16 


13 


Test  la 

Test  Ha 

Test  lb 

Test  Ic 

Test  lib 

Temperature  of  the  bottom 

1123  K 

1123  K 

1123  K 

Temperature  of  the  top  wall 

923  K 

1123  K 

923  K 

Temperature  of  the  side  walls 

Heat  flux  from  the  SOFC  stack 

773  K 

6  kW  m-2 

1123  K 

773  K 

lOkWm-2 

Thermal  conductivity  of  the  side  walls 

lWm-1K-1 

0.5  Wm-1  K-1 

‘Ambient’  temperature 

500  K 

500  K 

Fig.  2.  The  maximum  value  of  temperature  inside  the  reformer  as  a  func¬ 
tion  of  the  distance  between  the  fuel  reformer  and  SOFC  stack  and  side 
walls  emissivity,  for  the  first  problem  and  thermal  boundary  Test  la:  7b  (x, 
y,  0)=1123K,  7Y(x,  y,  zn)  =  923  K,  7s=773K,  thermal  conductivity  of  the 
reformer,  A,r  =  l.OWm-1  K-1  and  the  emissivity  of  the  horizontal  walls  was 
£  =  0.5. 


The  examples  of  numerical  results  for  the  first  problem  are 
shown  in  Figs.  2  and  3.  Fig.  2  shows  examples  of  the  results  of 
numerical  computations  for  thermal  boundary  Test  la.  The  tem¬ 
peratures  of  the  bottom,  top  and  side  walls  were  assumed  to  be 


Fig.  3.  The  maximum  value  of  temperature  inside  the  reformer  as  a  function 
of  the  distance  between  the  fuel  reformer  and  SOFC  stack  and  side  wall  emis¬ 
sivity,  for  the  first  problem  and  thermal  boundary  Test  Ha:  q  =  6000.0  Wm-2, 
To  =  500  K,  thermal  conductivity  of  the  reformer,  A,r  =  0.5  W  m-1  K-1  thermal 
conductivity  of  the  side  walls  =  l.OWm-1  K-1  and  the  emissivity  of  the 
horizontal  walls  was  s  =  0.5. 


constant;  the  temperature  of  the  bottom  is  7b  (x,  y,  0)  =  1 123  K, 
the  temperature  of  the  top  wall  is  7t(x,  y,  zn)  =  923  K,  and  the 
temperature  of  the  side  walls  is  7s=773K.  Fig.  2  presents 
the  maximum  temperature  inside  the  reformer  as  a  function 
of  the  distance  between  the  fuel  reformer  and  SOFC  stack 
and  the  side  wall  emissivity.  The  emissivity  of  the  horizontal 
walls  was  e  =  0.5,  the  thermal  conductivity  of  the  reformer  was 
Ar  =  l.OWm-1  K-1.  Fig.  3  shows  examples  of  the  results  of 
numerical  computations  for  thermal  boundary  Test  Ha.  The  heat 
flux  from  the  SOFC  stack  has  the  value  q  =  6000.0  W  m-2,  and 
‘ambient’  temperature  was  To  =  500  K;  the  thermal  conductivity 
of  the  side  walls  As  equals  1.0  W  m-1  K-1  and  thermal  conduc¬ 
tivity  of  the  reformer  Ar  equals  0.5  W  m-1  K-1 .  Fig.  3  presents 
the  maximum  value  of  temperature  inside  the  reformer  as  a  func¬ 
tion  of  the  distance  between  the  fuel  reformer  and  SOFC  stack 
and  side  walls  emissivity.  As  is  evident  in  Figs.  2  and  3,  the 
temperature  inside  the  reformer  is  strongly  affected  by  both  the 
position  of  the  module  chamber  and  by  the  walls’  emissivity.  As 
the  value  of  h  decreases  (i.e.  the  distance  between  the  reformer 
and  SOFC  stack  decreases)  the  temperature  inside  the  reformer 
increases.  The  calculated  temperature  distributions  inside  the 
reformer  show  that  the  value  of  the  temperature  depends  on 
the  walls’  emissivity  (i.e.  the  emissivity  of  the  bottom  reformer 
wall,  emissivity  of  the  SOFC  stack  and  the  side  walls  in  the 
SOFC  module).  Analyses  of  numerical  results  have  proved  to 
be  interesting.  For  example,  when  the  emissivity  of  the  side 
walls  decreases  the  temperature  inside  the  reformer  increases. 
However,  when  the  emissivity  value  of  the  bottom  reformer 
wall  and/or  the  emissivity  of  the  SOFC  stack  increase  then  the 
temperature  inside  the  reformer  increases  as  well  [26] . 

The  examples  of  numerical  results  for  the  second  problem  are 
shown  in  Figs.  4-9.  Figs.  4  and  5  show  examples  of  the  results 
of  numerical  computations  for  thermal  boundary  Test  lb.  The 
temperature  of  the  bottom,  top  and  side  walls  was  assumed  to  be 
constant;  7b (x,  y,  0)  =  7t(x,  y,  zn)  =  7s  =  1 123  K.  The  inlet  molar 
flow  rate  of  methane  has  the  value  ncu4  =  0.2  x  10  4  kmol  s  1 , 
the  steam- to-carbon  ratio  equals  2.5,  pressure  inside  the  reformer 
chamber  was  P  =  1.01  x  105  Pa,  the  catalysts  active  area  has  the 
value  Acat  =  10m2,  the  volume  of  reformer  equals  0.1  m3,  the 
distance  between  stack  and  fuel  reformer  was  h  =  0.35  m,  and 
the  emissivity  of  all  surfaces  in  the  system  £  equals  0.5.  Fig.  4 
presents  the  distribution  of  temperature  and  the  distribution  of 
methane  inside  the  reformer  in  the  longitudinal  direction  of  the 
fuel  reformer.  The  average  temperature  and  change  of  the  molar 
flow  rate  of  methane  in  the  longitudinal  direction  of  the  fuel 
reformer  are  illustrated  by  Fig.  5.  As  is  evident  in  Figs.  4  and  5, 
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Fig.  4.  Temperature  distribution  (A)  and  CH4  distribution  (B)  in  the  longitudinal  direction  of  the  fuel  reformer  (thermal  boundary  Test  lb). 
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Fig.  5.  Average  temperature  and  change  of  methane  molar  flow  rate  in  the 
longitudinal  direction  of  the  fuel  reformer  (thermal  boundary  Test  lb). 


Fig.  7.  Average  temperature  and  change  of  methane  molar  flow  rate  in  the 
longitudinal  direction  of  the  fuel  reformer  (thermal  boundary  Test  Ic). 


the  distribution  of  temperature  inside  the  reformer  is  not  uni¬ 
form,  and  when  the  molar  flow  rate  of  methane  decreases  the 
temperature  inside  the  reformer  increases  (Figs.  4  and  5).  The 
reforming  reaction  is  strongly  endothermic,  therefore  the  tem¬ 
perature  decreases  at  the  reformer  inlet. 

Figs.  6  and  7  show  the  results  of  numerical  computations  for 
different  configurations  of  thermal  boundary  Test  Ic.  The  tem¬ 
perature  of  the  bottom,  top  and  side  walls  was  assumed  to  be 


constant;  the  temperature  of  the  bottom  is  7b  (x,  y,  0)  =  1 123  K, 
the  temperature  of  top  wall  is  Tj(x,  y,  zn)  =  973  K,  and  the  tem¬ 
perature  of  the  side  walls  is  7s  =773  K.  The  inlet  molar  flow 
rate  of  methane  has  the  value  ^ch4  =0.1  x  10-4  kmol  s_1,  the 
steam- to-carbon  ratio  equals  2.5,  pressure  inside  reformer  cham¬ 
ber  was  P  =  1 .01  x  105  Pa,  the  catalysts  active  area  has  the  value 
Acat  =  10  m2,  the  reformer  volume  equals  0.1  m-3,  the  distance 
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Fig.  6.  Temperature  distribution  (A)  and  CH4  distribution  (B)  in  the  longitudinal  direction  of  the  fuel  reformer  (thermal  boundary  Test  Ic). 
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Fig.  8.  Temperature  distribution  (A)  and  CH4  distribution  (B)  in  the  longitudinal  direction  of  the  fuel  reformer  (thermal  boundary  Test  lib). 


between  stack  and  fuel  reformer  was  h-  0.1m,  and  the  emissiv- 
ity  of  all  surface  in  the  system  is  s  equals  0.5.  Fig.  6  presents  the 
distribution  of  temperature  and  the  distribution  of  methane  inside 
the  reformer  in  the  longitudinal  direction  of  the  fuel  reformer. 
The  average  temperature  and  change  of  the  molar  flow  rate  of 
methane  in  the  longitudinal  direction  of  the  fuel  reformer  are 
presented  in  Fig.  7.  As  is  evident  in  Figs.  4-7,  the  tempera¬ 
ture  inside  the  reformer  is  strongly  affected  by  the  position  of 
the  module  chamber  and  temperature  distribution  of  side  walls. 
As  the  temperatures  of  the  side  walls  decrease,  so  too  does  the 
average  temperature  inside  the  reformer. 

Figs.  8  and  9  exemplify  the  results  of  numerical  computa¬ 
tions  for  thermal  boundary  Test  lib.  The  heat  flux  from  the 
SOFC  stack  has  the  value  q  =  10,000  W  m-2  and  ‘ambient’  tem¬ 
perature  was  To  =  500  K;  the  thermal  conductivity  of  the  side 
walls  A,s  equals  0.5  Wm-1  K-1  and  the  emissivity  of  all  sur¬ 
faces  in  the  system  was  £  =  0.5.  The  inlet  molar  flow  rate  of 
methane  has  the  value  ncu4  =0.1  x  10-4  kmol  s-1,  the  steam- 
to-carbon  ratio  equals  2.5,  the  distance  between  stack  and  fuel 
reformer  was  h-  0.1  m,  the  pressure  inside  reformer  chamber 
was  P—  1.01  x  105  Pa,  the  catalysts  active  area  has  the  value 
Acat=10m2,  the  volume  of  reformer  equals  0.1m-3.  Fig.  8 
presents  the  distribution  of  temperature  and  the  distribution  of 
methane  inside  the  reformer  in  the  longitudinal  direction  of  fuel 


reformer.  The  average  temperature  and  change  of  the  molar  flow 
rate  of  methane  in  the  longitudinal  direction  of  fuel  reformer  are 
illustrated  by  Fig.  9.  Analyses  of  the  numerical  results  proved 
interesting.  For  example,  the  distribution  of  temperature  inside 
the  reformer  is  strongly  effected  by  the  kind  of  thermal  boundary 
adapted  in  numerical  tests. 

In  the  most  existing  literature  for  SOFC  model  the  radia¬ 
tion  is  often  simply  described  by  the  two  face-to-face  planes. 
An  example  of  comparison  of  numerical  results  obtained  by 
model  presented  in  this  paper  and  simplified  model  ‘two  face- 
to-face  planes’  is  shown  in  Fig.  10.  In  the  numerical  calculations 
presented  in  Fig.  10,  the  heat  flux  from  the  SOFC  stack  has 
the  value  q-  6000.0  Wm-2,  and  ‘ambient’  temperature  was 
To  =  500  K;  the  thermal  conductivity  of  the  top  wall  At  equals 
1.0 Wm-1  K-1  and  thermal  conductivity  of  the  reformer  Ar 
equals  2.0  Wm-1  K-1,  side  walls  were  insulated.  As  is  evident 
in  Fig.  10  the  calculated  temperature  distribution  depends  on 
the  numerical  modelling.  The  temperature  distribution  obtained 
by  simplified  model  is  much  lower  and  is  not  ‘realistic’  (i.e. 
we  cannot  see  maximum  value  of  temperature  in  central  part 
of  the  reformer).  Analyses  of  the  numerical  results  obtained 
by  these  two  models  have  proved  the  necessity  of  the  complex 
calculations. 


Fig.  9.  Average  temperature  and  change  of  methane  molar  flow  rate  in  the 
longitudinal  direction  of  the  fuel  reformer  (thermal  boundary  Test  lib). 


Fig.  10.  Temperature  distribution  at  the  bottom  wall  of  reformer  (comparison 
of  numerical  results  obtained  by  present  model  and  by  the  face-to-face  planes 
model). 
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4.  Conclusions 

The  present  paper  describes  the  numerical  modelling  of 
the  radiative  heat  transfer  process  in  the  module  chamber  of 
an  internal  indirect  reforming-type  SOFC.  As  in  any  high- 
temperature  system,  radiative  heat  transfer  effects  are  significant 
and  need  to  be  accounted  for  to  accurately  predict  tempera¬ 
tures  in  the  SOFC  module.  In  this  article,  governing  equations 
for  radiative  heat  transfer  are  described  using  Hottel’s  zone 
method.  The  steam  reforming  reaction  occurring  inside  the 
fuel  reformer  is  described  using  Achenbach  model.  The  result¬ 
ing  radiation-conduction-convection  conjugate  heat  transfer 
problems  are  numerically  solved  with  a  combination  of  the 
Gauss-Seidel  and  Newton-Raphson  methods.  The  obtained 
results  indicate  that,  for  the  development  of  effective  indirect 
internal  reforming,  the  position  of  the  reformer  in  the  module 
chamber  and  emissivity  of  the  surfaces  of  the  reformer,  cell  and 
other  elements  in  the  SOFC  module  all  play  a  key  role.  As 
the  distance  between  the  reformer  and  SOFC  stack  decreases 
the  temperature  inside  the  reformer  increases.  When  the  emis¬ 
sivity  of  the  side  walls  decreases  the  temperature  inside  the 
reformer  increases.  However,  when  the  emissivity  value  of  the 
bottom  reformer  wall  and/or  the  emissivity  of  the  SOFC  stack 
increases  then  the  temperature  inside  the  reformer  increases  as 
well. 
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